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Colloidal particles dispersed in a liquid crystal lead to distortions of the director field. The distortions are 
responsible for long-range effective colloidal interactions whose asymptotic behaviour is well understood. The 
short distance behaviour of the interaction, however, is sensitive to the structure and dynamics of the topolog- 
ical defects nucleated near the colloidal particles in the strong anchoring regime. The full non-linear theory 
is required in order to determine the interaction at short separations. Spherical colloidal particles with suffi- 
ciently strong planar degenerate anchoring nucleate a pair of antipodal surface topological defects, known as 
boojums. We use the Landau-de Gennes formalism in order to resolve the mesoscopic structure of the boojum 
cores and to determine the pairwise colloidal interaction. We compare the results in three {2>D) and two (2D) 
spatial dimensions. The corresponding free energy functional are minimized numerically using finite elements 
with adaptive meshes. Boojums are always point-like in 2D, but acquire a rather complex structure in 3D which 
depends on the combination of the anchoring potential, the radius of the colloid, the temperature and the LC 
elastic anisotropy. We identify three types of defect cores in 3D which we call single, double and split core 
boojums, and investigate the associated structural transitions. The split core structure is favoured by low tem- 
peratures, strong anchoring and small twist to splay or bend ratios. For sufficiently strong anchoring potentials 
characterised by a well-defined uniaxial minimum, the split core boojums are the only stable configuration. In 
the presence of two colloidal particles there are substantial re-arrangements of the defects at short distances, 
both in 3D and 2D. These re-arrangements lead to qualitative changes in the force-distance profile when com- 
pared to the asymptotic quadmpole-quadrupole interaction. In line with the experimental results, the presence 
of the defects prevents coalescence of the colloidal particles in 2D, but not in 3D systems. 

PACS numbers: 61.30.Dk,61.30.Jf,83.80.Xz 



1. INTRODUCTION 

Controlling the self-assembly of colloidal particles is an im- 
portant aspect of pure and applied colloid science. In this 
context, assembling novel structures in liquid crystal matri- 
ces has become a very active field of research as a result of 
important theoretical and experimental advances Ill-Ill. Liq- 
uid crystals (LCs) are characterized by anisotropic mechani- 
cal and optical properties due to the long-range orientational 
molecular ordering. Consequently, in addition to the usual 
isotropic colloidal interactions, colloids in liquid crystals in- 
teract through strongly anisotropic effective forces, and have 
been observed to self-assemble into ordered aggregates in two 
(2D) 10] and three (3D) spatial dimensions 15]. A distinctive 
feature of these systems is the presence of topological defects, 
the type of which determines the symmetry of the far-field di- 
rector around an isolated colloidal particle, which in turn gov- 
erns the asyrnptotic behavior of the nematic-induced colloidal 
interactions |l6i]. The behavior of the interactions at short dis- 
tances is more complex since the defects can move when the 
colloidal particles are brought close together, rendering the 
description of the self-assembly of liquid crystal colloids a 
challenging theoretical problem. 
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Theoretical and simulation methods have been developed 
to study LCs on various scales. On the macroscopic scale, 
the long-range forces between colloidal particles can be calcu- 
lated using the Frank-Oseen (FO) elastic free energy approach 
||6t]. On the mesoscopic scale where non-Hnear effects be- 
come important, a description based on the Landau-de Gennes 
(LdG) tensorial order-parameter theory is necessary [Q. This 
is the case when the interactions between the defects domi- 
nate, as revealed by the changes in their equilibrium structure 
when the distance between the paiticles decreases [8]. On the 
microscopic scale, presmectic ordering and density variations 
cannot be ignored and density functional approaches Q or 
computer simulations are usually employed lUOll . 

The anchoring at the surfaces of colloidal particles can be 
controlled, for instance, through the adsorption of surfactants 
or deposition of thin organic or inorganic films lITTl - ITil . When 
the anchoring on a spherical colloid is planar, a tangential ne- 
matic director configuration is expected at the colloidal sur- 
face. At the FO level the boundary conditions are met by 
the creation of a pair of antipodal surface topological defects, 
called boojums ifTsl [l6ll . The boojums are aligned with the 
far-field director, ensuring that the far-field LC configuration 
is homogeneous as required by the global uniform boundary 
conditions UtIi . 

Recently, experimental measurements of the effective 
forces between two spherical colloids with planar anchoring 
in a 3D nematic have been reported |[l3l . The authors have 
found that the forces at short distances deviate considerably 
from the asymptotic quadrupole-quadrupole ones, and that the 
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equilibrium configuration of a pair of particles corresponds to 
close contact implying the absence of a repulsion at short dis- 
tances. By contrast, experimental results obtained for inclu- 
sions in free-standing smectic C films ifTsl [l9ll revealed the 
presence of a short distance repulsion which keeps the inclu- 
sions at a well-defined equilibrium separation larger than con- 
tact. Again, as in the 3D system of Ref. 1IT3I1 . the asymptotic 
quadrupole-quadrupole approximation for the pair-wise force 
breaks down at short distances. 

A recent theoretical study based on the numerical mini- 
mization of the 3D LdG free energy functional 1I20I1 is in line 
with the experimental results of Ref. ifisll concerning the equi- 
librium configuration of two colloidal spheres and the absence 
of colloidal short- rang e repulsion. Theoretical analysis of 
the 2D FO model llSl fl9ll . based on an ad hoc ansatz for a 
single particle solution and the superposition approximation, 
have correctly described the configuration of two inclusions 
in a smectic C film at short distances, but have failed to pro- 
vide a consistent mechanism for the short-range repulsion ob- 
served experimentally. A subsequent superposition analysis 
ll2lll . based on the exact single particle solution ll22ll . has pre- 
dicted the existence of two regimes as a function of the an- 
choring strength. For weak anchoring the inclusions coalesce 
while for strong anchoring, the effective interaction exhibits a 
well defined minimum and a strong short-distance repulsion 
that prevents coalescence 1I21I1 . 

Although the 3D studies of Ref. ll20ll are in qualitative 
agreement with the experimental results of Ref. lITsll . the 
question of the structure of the boojum cores remains open. 
The authors of Ref. ll20ll could not resolve the structure the 
defect cores because the meshes which have been used in 
the minimization of the LdG free energy functional were too 
coarse. In the 2D systems, studied in Refs. ifisi fT9l l2lll . the 
theoretical results are far from satisfactory since the superpo- 
sition approximation is not valid at short distances where the 
interaction between defects is expected. 

Motivated by these open questions, and by the experimen- 
tal results which reveal a quaUtative difference between 2D 
and 3D systems IH [H [H HI, we address in this article 
the issues of the structure of the boojum cores and the ef- 
fective interactions between circular and spherical colloids in 
2D and 3D, respectively. We concentrate on the mesoscale 
and present the results of systematic numerical analysis of 
the LdG free energy functionals in 2D and 3D for colloidal 
particles characterized by a tangential (degenerate in 3D) an- 
choring. We use finite element methods (FEM) with adaptive 
meshes in order to minimize the LdG functionals and obtain 
pair-wise effective interactions at all distances. We find that 
the defects are point like in 2D but acquire a rather complex 
structure in 3D. The structure depends on the combination of 
the anchoring potential, the particle radius, the temperature, 
and the LC elastic anisotropy. More specifically, we identify 
defect core transitions between, as we call them, single-core, 
double-core, and split-core structures. The split-core structure 
is favoured by increasing the anchoring strength or the parti- 
cle radius and by decreasing the temperature. The LC elastic 
anisotropy responsible for a preferred planar anchoring at the 
nematic-isotropic interface (twist elastic constant smaller than 



bend and splay) ||24|l also favours the split-core structure. Fi- 
nally, the split-core structure is favored by surface potentials 
characterized by a well-defined uniaxial minimum llzsll . In 
2D the boojums have a single-core point-like structure for ar- 
bitrarily strong surface potentials. 

The paper is organized as follows: In Sec. |2]we discuss the 
3D system. First we review the experimental results and de- 
fine the LdG functional. Then we discuss the results for a sin- 
gle spherical colloid, focusing on the structures of the boojum 
core as a function of the anchoring strength, the colloid size 
and the temperature. We consider two types of LCs differing 
by their anisotropy. Finally we consider the effective interac- 
tion between two colloidal particles and compare the results 
with the experimental measurements |[l3l I23I1 . In Sec. |3]we 
turn our attention to 2D systems and calculate the interaction 
between colloids, where a short-range repulsion is found for 
strong anchoring. A comparison with the experimental results 
of iflsi [19I1 is carried out. In Sec. |4]we discuss and compare 
our results in 2D and 3D. For completeness, the numerical 
techniques used throughout this article are described in Ap- 
pendix [A] 



2. THREE-DIMENSIONAL SYSTEMS 

In the linear regime the effective force between two col- 
loidal particles can be computed within the one-elastic- 
constant approximation by using the electrostatic analogy 
For the case of spherical colloids with degenerate tangential 
anchoring the leading asymptotic term in the multipole ansatz 
for the director field is given by a quadrupole term r^^ |2^, 
where r is the distance to the center of the colloidal particle. 
The superposition approximation then leads to the quadrupo- 
lar effective interaction between two particles [ Ij]. 

Smalyukh et al. ifisll have studied the anisotropic interac- 
tions between two colloids with tangential anchoring by using 
laser tweezers. At (relatively) large distances d the radial de- 
pendence of the measured force was found to comply with 
the expected d^^ quadrupolar behaviour However, the angu- 
lar dependence of the force disagreed qualitatively with the 
quadrupolar one at all distances, indicating the importance of 
the non-linear effects and a breakdown of the superposition 
approximation. 

In what follows, we shall use the LdG theory in order to 
calculate effective interactions between two colloids with tan- 
gential anchoring. We shall consider two types of LCs dif- 
fering by their elastic anisotropics. Strong emphasis is placed 
on short distance behaviour and strong anchoring limit, where 
the non-linear effects are expected to dominate. The prob- 
lem is challenging from the theoretical point of view because 
the defect structure and dynamics are complex, both in the 
single-particle case as well as for interacting particles, as it 
was predicted a decade ago HI in 2D and described recently 
in 3D ll20ll . In the following we shall show that the colloidal 
interaction at short distances results, not only from the re- 
arrangement of the defect positions, but also from structural 
changes of the defect cores. 

Ultimately, we aim at describing quantitatively the experi- 



3 



mental results of Smalyukh et al. ifisll . for the angular depen- 
dence of the elastic force at moderate and short inter-particle 
distances. 



2.1. Landau-de Gennes free energy functional 

Within the Landau-de Gennes (LdG) theory iIztIi nematic 
liquid crystals are characterised by a traceless symmetric 
order-parameter tensor with components Qij, which can be 
written as 

Qij = {run j - ^ dij^ + {Ulj - minij) , (1) 

where n,- are the Cartesian components of the director field n, 
Q is the uniaxial order-parameter, which measures the degree 
of orientational order along the nematic director, and B is the 
biaxial order parameter, which measures the degree of orien- 
tational order along the directions perpendicular to n, charac- 
terized by the eigenvectors 1 and m. The corresponding LdG 
free energy functional is 

^LdG= / {fh+fel)d\x+ f f,ds (2) 
Ja JdQ. 

with fi, and /g/ the bulk and elastic free energy densities, given 
by 

ft, = aTr(f-bTr(f + c{Tr(f)\ (3) 
/./ - y'^A-e/y-'JAGo + y'^iGu^^Gi-t: (4) 

where a depends linearly on the temperature T and is usu- 
aly written as a = aQ{T — T*), with ao a material dependent 
constant and T* the supercooling temperature of the isotropic 
phase, b and c are positive (material dependent) constants, 
and Li and L2 are phenomenological parameters which can be 
related to the Frank-Oseen (FO) elastic constants. The first in- 
tegral in Eq. (|2]i is taken over the 3D domain, i2, occupied by 
nematic, while the second integral is over the surfaces dQ. (in 
our case the surfaces of the colloidal particles) and accounts 
for non-rigid anchoring conditions. 

Depending on the preferred orientation of the director (an- 
choring direction) with respect to the surface normal v, the 
surface anchoring is i) homeotropic, when the anchoring di- 
rection is parallel to v, ii) planar, when the anchoring direc- 
tion is orthogonal to v, and iii) tilted, when the anchoring 
direction and v form an angle smaller than n/2 (and larger 
than zero). The last two cases can be classified futher as 
monostable, multistable or degenerate (in some articles the 
term "random" is used instead of "degenerate"), depending on 
whether the surface imposes one, a finite or an infinite num- 
ber of equivalent anchoring directions lfl2ll . respectively. The 
simplest quadratic surface free energy, favouring monos- 
table nematic ordering Q^ i.e., with a well-defined director, 
scalar and biaxial order-parameters was proposed by Nobili 
and DurandlHi]: 

/. = w(G,/-G/,)' (5) 



When W > 0, the surface free energy has a unique mini- 
mum at Q = Q'' . In general this is not necessary, the require- 
ment being that the total free energy is bounded from below. 
In this paper we consider only planar degenerate anchoring 
described by a family of covariant surface potentials originally 
proposed by Fournier and Galatola llzsl : 

/v = Wi (g„ -Qijf + W2 (qI - (^) ') (6) 

where Q,j = Qij + G/,%, Qfj = {8,k - VfVk) Qij {5,j - ViVj), 
W\ is the anchoring strength favouring tangential orientation 
of the director n, and W2 > guarantees the existence of a 
minimum for a scalar order-parameter at the surface equal to 
its bulk value, Gi- At a flat surface the nematic is uniform 
and uniaxial everywhere. In the original formulation ll25ll the 
surface scalar order-parameter is allowed to vary and the bi- 
axiality at a flat surface increases with the difference between 
the surface and bulk scalar order-parameters. The quartic sur- 
face potential given by Eq. ^ may be viewed as the minimal 
biaxiality potential characterised by a well-defined degener- 
ate planar minimum. The quadratic surface potential (W2 — 0) 
is the covariant version of Eq. (|5]l the minimum of which de- 
pends on the coupling to the bulk nematic. 

In the strong anchoring regime, the nematic director is par- 
allel to the spherical colloidal surface everywhere, with an 
orientation determined by the far field director In Sec. 12.21 
we will see that the quartic term W2 has a profound effect on 
the structure of the topological defects, by controlling the cou- 
pling to the bulk nematic through the deviation of the surface 
and bulk scalar order-parameters. In the flat surface Umit the 
biaxiality vanishes in the strong anchoring regime [l25l|, while 
on the surface of large spherical colloids point-like singular- 
ities with charge + \ split into pairs of point-like singularities 
with charge H-1/2 connected by disclination lines. 

It is convenient to define the dimensionless temperature 
T = 24flc/fe^. At T < 1 the uniaxial nematic is stable and the 
degree of orientational order is given by 

The nematic becomes unstable atT>9/8. AtT=l both 
the nematic and the isotropic phases coexist. Typical values 
of the bulk parameters for 5CB are S ao = 0.044 x 10^ 
J /Km^ b = 0.816 X 10^ J/m^, and c = 0.45 x 10*' J/m^ 
Li = 6 X \Q-^^]/m, T* = 307 K. The spatial extension of 
inhomogeneous regions and the cores of topological defects 
is of the order of the bulk correlation length, which is given 

by ^ = (8c(3Li +2L2) /b^)^^^ at the nematic-isotropic (NI) 
transition ll30ll . 

We define the elastic constant anisotropy rj = L2/L1, and 
consider two cases i) 77 = 2 with £, ~ 15nm corresponding to 
5CB, and T] = — 1/2 with £ ~ 8nm. Note that stability argu- 
ments require 77 > —3/2 IBlll . Depending on the LC material, 
L2 can be positive or negative, and its sign controls the molec- 
ular orientation at the NI interface, which is planar for T] > 
and homeotropic for 77 < ll24ll . 
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The LdG elastic constants Li and L2 may be related to the 
FO elastic constants 13211 . Ki = Ko, and K2, through the uniax- 
ial ansatz Qij = (3/2) Qh {tiiiij — 5,7/3), yielding Ki = ^ 
9Ql{Li +L2/2)/2 and K2 = 9QILi/2. In general Ki and K3 
are different, but in most cases the difference is small and the 
LdG free energy is deemed adequate. 

In the following the LdG free energy Eq. (|2]i is minimised 
numerically using FEM with the adaptive mesh refinement. A 
detailed description of the numerical procedures is given in 
AppendixlAl 



2.2. Single spherical particle: three types of boojum cores 

When spherical colloids with strong planar anchoring are 
dispersed in nematic LCs a pair of topological surface defects, 
named boojums, appears at the antipodes of the particles ITsl - 
[T7I1 . The understanding of the structure of defects in terms of 
the order-parameter distribution in the cores is beyond mere 
topological arguments and has been the subject of investiga- 
tion since the early 1930s lf33ll . For infinitely strong anchor- 
ing and within the class of axially symmetric fields, hedgehog 
defects exhibit three different structures: the radial hedge- 
hog, a small ring or loop disclination ||34ll and a third struc- 
ture, a split-core defect, which was found to be metastable 
1I35I1 . Phase and bifurcation diagrams indicate that the tran- 
sition from the hedgehog to the ring structure is first order, 
as predicted in Ref. [34]. Although several studies addressed 
the structure of hedgehog defect cores ||34|-|36|| the structure 
of boojums has not been fully understood OTIl . 

A homotopy classification of surface topological defects is 
given by Volovik in Ref. fJSi). According to Volovik's classi- 
fication a point defect at the boundary of a nematic LC may 
be viewed as the combination of a bulk hedgehog and a sur- 
face boojum and is therefore characterised by two topological 
charges - the charge of the bulk hedgehog and the index m 
of the projection t = n — v(n • v) of the director field onto the 
surface with normal v. The charges and m are related to a 
continuous topological charge .e/ defined as Ii39t i40l1 
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where a is a hemisphere surrounding the defect, and Q , (j) 
are arbitrary coordinates on a. Evaluation of the integral 
in Eq. ^ gives an explicit relation between the continuous 
charge ^ and the charges and m 1331 



m , 
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which provides a means of calculating the charge as the in- 
dex m can be calculated independently ll40ll . If si is an integer 
the defects can detach from the surface and if ^ = they may 
vanish altogether By contrast, for non-integer charges si the 
defects are "topologically" bound to the surface. For instance, 
a single-core boojum, which is discussed below, has si — 
and N i.e., the defect cannot vanish nor detach itself from 
the surface. 



In this section we investigate the dependence of the boojum 
core structure on i) the strengths of the anchoring potential, 
WijVVi, ii) the reduced temperature T, and iii) the colloidal 
radius R. We introduce the dimensionless anchoring strengths 
Wi = WiOlR (' = 1 : 2), and consider LCs with positive (e.g. 
5CB) and negative elastic anisotropies. We pay particular at- 
tention to the nature and degree of the nernatic order within 
the boojum cores. In agreement with Ref. [ 41: 1 we find that the 
apparent singularities of the director field are replaced by uni- 
axial order-parameter distributions with negative scalar order- 
parameter, corresponding to oblate nematic order, surrounded 
by biaxial layers. 

We shall show that boojum cores are axially symmetric 
point-like with index m = +1 on small colloids. We call this 
structure single core boojum. It is stable at high temperatures 
and relatively weak anchoring. In Ref. jilj a similar boo- 
jum core structure (named by the authors "fingered" boojum) 
has been reported at a flat surface. On large colloids, at low 
temperatures and strong anchoring, the axial symmetry is bro- 
ken and the boojum +1 point-like cores split into pairs of +5 
point-like surface defects which are connected by disclination 
lines. A structure without a fully developed disclination line, 
the double core boojum, is also found. The detailed structure 
of the boojums as well as the transitions between the differ- 
ent configurations depend in detail on the colloid and the LC 
parameters as we shall discuss below. 

We start by describing the three distinct stable configura- 
tions of boojum cores, namely single, double and split core 
boojums. The cores differ both in their surface and bulk struc- 
tures. Typical core configurations are illustrated in the lower 
panels of Fig. [1] where the degree of biaxiality defined as 
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is shown by a color map. )3 = characterises the uniaxial ne- 
matic while j8 = 1 corresponds to the maximal biaxiality ne- 
matic state. j3 = 1 is obtained when one (and only one) eigen- 
value of Q vanishes. Note that in the isotropic phase all of the 
eigenvalues of Q vanish Issll . Figure [T]illustrates the effect of 
T on the structure of boojums for rj =2 and strong quadratic 
(w2 — 0) surface potential. As T decreases the boojum core 
transforms from the single to the double core, and then to the 
split core structure with the accompanying disclination half- 
ring. The core is fairly large close to the NI transition and 
shrinks as the temperature decreases. 

A description of the boojum core configurations in terms of 
the eigenvalues of Q along the z-axis (the direction of the far- 
field director) is illustrated in Fig. |2] where the eigenvalues 
~ 1,2,3 are plotted as functions of the distance z to the 
colloidal surface at fixed x = y = (the x,y coordinates of the 
center of the colloidal particle). 

Single core boojums are uniaxial along the z— axis, as re- 
vealed by the two degenerate eigenvalues. The scalar order- 
parameter vanishes at z/^ ~ 0.5 above the colloidal surface, 
where the three eigenvalues vanish and the system is isotropic. 
Below this point the scalar order-parameter is negative, since 
the principal (non-degenerate) eigenvalue is negative and the 
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a b c 

FIG. 1: Three types of boojum core structures: a) single core, b) dou- 
ble core, and c) split core configurations. The transitions between 
different configurations are driven by decreasing temperature t: a) 
T ~ —0.3, b) 1 ~ —0.7, c) T ~ —3.2, at fixed anchoring strengths, 
w\ ~ 37,H'2 = 0, radius of the colloidal particle, R = 0.1/im ~ 6.6i^, 
and elastic constant anisotropy, rj = 2. The scale bars represent 
the nematic correlation length ^ at the NX transition. In the upper 
panels and in the lower panel c) the iso-surface of the scalar order- 
parameter corresponding to Q = 03Qh is shown. In the lower panel 
b) the proto-ring of the double core boojum is represented by the 
iso-surfaces Q = OAQj, and Q = 0.22/,. The lower panels depict the 
color-coded biaxiality parameter Eq. ( llOt in the plane z = R- The 
biaxiality parameter in the plane y = is also shown for the single 
core boojum, in the lower panel a). The short lines on the surface of 
the colloidal particle represent the eigenvector of Q corresponding to 
the largest eigenvalue. This eigenvector corresponds to the nematic 
director only in the case of prolate nematic order. See the main text 
for further details. 

nematic order is oblate. Above z/S, ~ 0.5 the scalar order- 
parameter is positive indicating that the nematic becomes pro- 
late. 

By contrast, double core boojums are biaxial along the 
z— axis. The surface of maximal biaxility intersects the z— axis 
at three distinct points (see Fig.lU? and Fig.O separating re- 
gions of low nematic oblate and prolate order. 

Split core boojums are also biaxial along the z— axis. The 
surface of maximal biaxiality intersects the z— axis at two dis- 
tinct points, z/^ — 0.49 and z/^ — 0.63, that delimit the re- 
gion of oblate nematic order. In the outer region, z/^ > 1, 
two of the eigenvalues are degenerate and the nematic liquid 
crystal is uniaxial. As z decreases the two negative eigenval- 
ues depart from each other, implying non-zero biaxiality. One 
eigenvalue remains negative. The other two exchange places 
at a point that can be identified with the center of the boojum 
core. At this point the nematic is uniaxial with Qi, < 0. All 
such points form a 3D curve which is surrounded by the sur- 
face of maximal biaxiality (detQ = 0). This surface encloses 
the region with oblate uniaxial order {Qb < 0). 

In all cases, the boojum cores are non-singular and biaxial 
over extended regions, with linear dimensions of the order of 
the bulk correlation length, ~ ^ . We proceed to describe in 
detail the structure of the three types of boojum cores. 

Single core boojums are characterised by an axially sym- 
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FIG. 2; Eigenvalues of the tensor order-parameter Qij on the z— axis 
for the three types of boojum cores: a) single core, b) double core, 
and c) split core configurations, respectively. The model parameters 
are as in Fig.[T] The far-field director is oriented along the z— axis. 




FIG. 3: Structure of the double core boojum represented by the iso- 
surface detQ = 0, and the color-coded biaxial parameter Eq. ( 110) 
in the plane z = R- a: the model parameters are as in Fig. QJ'; b: 
T~0.16,i?~9.3^ 



metric distribution of the order-parameter Q. Consequently, 
j8 = on the z axis, where two eigenvalues are degenerate, 
as shown in Fig|2]3. The iso-surface detQ = is a distorted 
hemisphere that is biaxial ()3 = 1) everywhere (See FigUJi) 
except at the point where it intersects the z axis, z/£, — 0.5, 
where the fluid is isotropic. At z = 0, the location of the pu- 
tative surface point defect, the scalar order-parameter is nega- 
tive, i.e., the molecules are perpendicular to the director 112711 
(oblate nematic order). 

In double core boojums the axial symmetry is broken. As 
a result the system becomes biaxial along the z axis, Figs.[Tj7 
and The most prominent feature of this structure how- 
ever, is the appearance of a proto-ring resembling the fully 
developed half -ring disclination of the split core boojum. The 
proto-ring is formed in a region, of the order of ^, of low 
nematic order (the scalar order-parameter is < 0.52/,) and is 
depicted in Fig.[Tl7 by the iso-surface Q ~ O.lQi,. The double 
core structure shares a number of features with the split core 
boojum. Its innermost core is a half-ring starting and end- 
ing on the colloidal surface, with very low biaxiality, j3 ~ 0, 
and negative scalar order-parameter. Further inspection of the 
eigenvalues of Q reveals that the iso-surface detQ = (j3 = 1) 
is a hemisphere with a handle at the top, see Fig. |3] This 
iso-surface intersects the colloidal surface in a single closed 
loop. In the hole, below the handle, the biaxiality is vanish- 
ingly small and the scalar order-parameter is positive, i.e., the 
nematic is prolate uniaxial. 

Finally, the split core boojum is characterised by a ^ discli- 
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nation half-ring connected to the colloidal surface. The inner 
core of the ring has zero biaxiality and oblate nematic order. 
This region is surrounded by the iso-surface det Q = 0, form- 
ing a half-torus see Figs.[T]: and|2]:, which intersects the col- 
loidal surface in two distinct closed loops. The nematic is 
prolate uniaxial everywhere except in a region of order ~ ^ in 
the core of the half-ring. 

How does the stability of these structures depend on the 
colloidal and LC parameters? All three types of boojums are 
topologically equivalent, in the sense that any of them can be 
transformed smoothly into any other The split and double 
core boojums result from breaking the axial symmetry of the 
single core, which splits on the colloidal surface. This split- 
ting is complete in the split core boojum in the sense that it 
is accompanied by the splitting of the detQ = surface, and 
only partial in the double core boojum. The dependence of 
the equilibrium (stable) boojum core structure on the strength 
of the anchoring potential, , W2, the reduced temperature T, 
and the radius R is complex and will be addressed below for 
LCs with positive and negative elastic anisotropies. 




2.2.1. Strength of the anchoring potential 

Usually the anchoring on the surface of the colloidal par- 
ticle is assumed to be infinitely strong. Here we investigate 
in detail the effects of the strength of the quadratic (wi) and 
quartic (W2) terms of the surface potential, Eq. (|6]l. When the 
anchoring strengths wi and W2 are large the split core structure 
is favoured independently of R,T, and 77. When W2 = and 
wi is large the split (double) core boojum is found only for 
sufficiently large colloids and low temperatures in LCs with 
T] > (77 < 0). At flat surfaces the auartic potential induces 
uniform uniaxial nematic profiles ll25ll . On spherical colloids, 
in the strong anchoring regime, this potential stabilizes the 
split core boojums where the biaxial regions are reduced. By 
contrast, in systems with a quadratic surface potential only, 
W2 = 0, the minimum is determined by the coupling of the 
surface to the bulk nematic Izsll . In this case the stability of 
the boojum cores depends in detail on Wi, R, T, and 77. 

We start by investigating the quadratic surface potential, 
W2 ~ 0. In Fig|4]we plot the elastic free energy F^'^ — Fg, 
where F^'^ = ^Ldo — fbiQh)^ is the excess (over the uniform 
uniaxial nematic) free energy and is the surface free en- 
ergy, as a function of the dimensionless anchoring strength, 
w = wi, for several values of T from just below the NI tran- 
sition (t ~ 0.97) to deep in the bulk nematic (t ~ —1.6), at 
fixed values of rj and R. Q. is the LC volume. 

Two regimes are observed. In the strong anchoring regime, 
w ^ 1, the elastic free energy exhibits a plateau that starts 
with the nucleation of the boojums at opposite poles of the 
colloidal particle. The inset shows that in this regime the ex- 
cess, F^^\ and the elastic, F''' — F^, free energies approach 
each other asymptotically, since the surface free energy de- 
creases roughly as ' . In the weak anchoring regime, w < 1 , 
the surface free energy increases linearly with w and domi- 
nates over the elastic term, which grows as w^. This behaviour 
is in line with the case of the homeotropic anchoring, where 



FIG. 4: Elastic free energy (excess over the uniform nematic free 
energy minus the surface free energy Fs) as a function of the 
dimensionless quadratic anchoring strength, w = wi, {wo =0), at 
several values of reduced temperature T. The inset shows the elastic 
(full line) excess (dashed line) and surface (dash-dotted line) free en- 
ergies as functions of w. The lower panel depicts the configurations 
of (1) the single core boojum for w ~ 3.7, and (2) the split core boo- 
jum for w ~ 37, T ~ — 1 .6. The black lines on the surface of the col- 
loidal particle represent the eigenvectors corresponding to the largest 
eigenvalue of Q. The iso-surfaces coiTesponding to Q = 0.3Qi, are 
also shown. The LC anisotropy Tj = 2, and the radius of the colloid 
R/^~6.6. 



the elastic free energy was also found to vary as in the 
weak and to saturate in the strong anchoring regimes ll42l . 

Figures Hfl) and (2) illustrate how the structure of the boo- 
jum cores evolves as the strength of quadratic surface poten- 
tial increases, at low temperatures. As the anchoring strength 
increases the competition between elastic and anchoring en- 
ergies leads to conflicting director orientations at the poles of 
the colloidal particles where the regions of reduced uniaxial 
and increased biaxial order develop. Inspection of the direc- 
tor configurations confirms that the system starts by develop- 
ing a single core point-like surface defect, see Fig|4jl). The 
projection t of the director field onto the surface of the par- 
ticle is characterised by an index m = +1. In this regime of 
small w, the surface director field has a finite normal compo- 
nent as can be seen in FigHJl). Deep in the strong anchoring 
regime the normal component of the surface director vanishes, 
see Fig|4l[2), and the region of reduced uniaxial nematic order 
opens into a half-ring which is attached to the colloidal sur- 
face, and resembles a j disclination line defect. 

For the quadratic surface potential (w2 = 0) boojum core 
structural transitions driven by the anchoring strength are ob- 
served at low temperatures only. The boojum core structure 
diagram is shown in Fig.|5]by circles. The single core is sta- 
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FIG. 5: Defect core transitions driven by the anchoring strength w. 
Circles correspond to a quadratic anchoring potential, wo = 0, while 
squares and diamonds correspond to a quartic anchoring potential, 
w\ and W2 7^ 0. The surface core splitting is the distance between 
two ^ point-like surface defects in the double and split core boo- 
jums. Single core boojums are stable at small anchoring strengths 
and high temperatures. Split core structures are stable at large an- 
choring strengths and low temperatures. 7] = 2, and the colloid ra- 
dius R/^ — 6.6. Solid symbols correspond to double core boojums 
that are stable over a narrow range of anchoring strengths. 



FIG. 6: Defect core transitions driven by the temperature. R = O.ljlm 
which is ~ 6.6^ for T] =2, and ~ \2A^ for 7] = — ^. Split core 
boojums are stable for the quartic, w\ and wo 7^ 0, potential at all 
temperatures and for both r\ =2 (diamonds) and T] = — (squares). 
For the quadratic potential, W2 = and LCs with rj =2 the sequence 
single-double-split cores is observed as the temperature decreases. 



Strengths, for both types of surface potential. 



ble at weak anchoring w < 6, the double core is stable in a 
narrow range of anchoring strengths, 6 < w < 15, while the 
split core is stable at strong anchoring w > 15. At high tem- 
peratures (t = 0.16,0.97, not shown in Fig|5]l the quadratic 
anchoring potential is incapable of stabilizing the double or 
spilt core boojums. Only the single core structure is observed 
in the range of anchoring strengths 0.00037 < w < 370. 

We now proceed with the analysis of the quartic potential, 
Eq. (|6]l. We set 77 = 2, /? ~ 6.6^ and assume, for simplic- 
ity w = wi = W2- We find three types of stable cores, at 
all temperatures. As the strength w of the anchoring poten- 
tial increases, the sequence of single, double and split core 
structures is observed. Two cases corresponding to T ~ 0.16 
(squares) and T ~ 0.97 (diamonds), are shown in the core 
structure diagram of Fig. |5] At high temperature the single 
core is stable for w < 4.9 while the split core is stable for 
w > 5.5. The double core is stable for intermediate anchor- 
ing strengths w. The transition between these core structures 
occurs at lower values of w and the surface core splitting de- 
creases as T decreases. When compared with the quadratic po- 
tential at the same temperature and colloidal radius, the tran- 
sitions for the quartic potential occur at lower w. 

In summary, at fixed colloidal size and temperature, a 
strong quartic potential stabilizes the split core structure. A 
strong quadratic potential is capable of stabilizing the split 
core boojums, only for LCs with positive elastic anisotropy. 
For LCs with negative elastic anisotropy the split core boo- 
jums do not appear. In this article we address only two cases 
rj =2,7] = — 7. Presumably, there exist a positive threshold 
value of 77 above which quadratic surface potential favours the 
split core boojums. The double core structure is an interme- 
diate structure that is stable over a limited range of anchoring 



2.2.2. Temperature 

In the previous section we have investigated the role of the 
anchoring strength on the stability of the boojum core struc- 
tures, and found evidence that the temperature plays a major 
role on the stability of the boojum cores. Here we consider the 
temperature as the control parameter of the structural booj um 
core transitions. 

We start by investigating a colloidal particle with a 
quadratic surface potential, W2 = 0. We consider reduced tem- 
peratures in the range —3.24 < T < 1 and set R = O.ljlm, 
which is ~ 6.6^ in LCs with positive elastic anisotropy rj =2, 
and ~ 12.4^ in LCs with rj ~ —j, since the bulk correlation 
length £, depends on the value of 77 . 

A quadratic potential (w2 = 0) with wi ~ 37 stabilizes the 
split core structure at low temperatures only in LCs with pos- 
itive elastic anisotropy, 7j = 2. The core structure diagram is 
depicted in Fig. |6]by circles. For T < — L12 the split core 
structure is stable, while at high temperatures T > —0.27 the 
boojum has a single core structure. The transition between 
these structures seems to occur continuously through the inter- 
mediate double core structure, which is observed in the range 
— L12 < T < —0.27. The surface splitting of the split core 
structure increases slowly as the temperature decreases and 
saturates at a value close to ^ . 

The split core structure was not observed in LCs with neg- 
ative anisotropy, 77 = — 5 . Recalling that 77 > corresponds 
to Ki=K^> K2, while 77 < to /Ti = /Ts < K2, we note that 
unfavourable twist distortions are responsible for the suppres- 
sion of the split core structures in LCs with 77 = —5. Twist 
deformations vanish when the azimuthal component of the di- 
rector field, n^, vanishes ll43ll . which is the case for single core 



boojums. By contrast, for double and split core boojums nip 
acquires a non-zero value that increases from the double to the 
split core structures. Because the twist distortions are energet- 
ically unfavourable in LCs with 77 < the split core boojums 
are unstable with respect to the double core structure. This 
is in line with a twist transition that occurs at the core of a 
hyperbolic hedgehog defect (see Figs. 4, 5 in Ref. ll43ll ') on 
a colloidal particle with homeotropic anchoring, and in gen- 
eral, with the effect of elastic anisotropy on the structure of 
distortions around topological defects in 2D H4\ . 

For a quartic potential with w\ =W2 — 37 the split core is 
stable for both types of LC elastic anisotropy, at all tempera- 
tures. The core splitting increases with temperature up to 2.8^ 
for T] = 2, and 4.9^ for 77 = —5, close to the NI coexistence 
temperature (see Fig. |6] diamonds and squares, respectively). 
The quartic term in the surface potential is quite effective in 
stabilizing the spUt core structure, even when twist distortions 
are unfarourable (tj = —5). In fact, the split core is not only 
stabilized, but its splitting is even larger than that observed 
in systems with 77 > 0. We note that the boojum size is pro- 
portional to the splitting parameter shown in Fig. |6]and it de- 
creases as the temperature decreases. 

In summary, strong quadratic surface potentials and low 
temperatures favour split core boojums only for LCs with pos- 
itive elastic anisotropy. When 77 < the double core struc- 
ture is stable. Strong quartic potentials favour split core boo- 
jums, which are the only stable structures, regardless of the 
LC anisotropy. 



2.2.3. Particle radius 

Finally, we investigate the effect of the colloidal radius R on 
the structural transitions of the boojum cores. For this purpose 
we consider colloids with sizes in the range Q.OSjim ^ R ^ 
O.Sjim. 

Boojum core transitions driven by the colloidal size, for 
quadratic (circles) and quartic (diamonds) potentials are illus- 
trated in Fig. Q For LCs with f] ~2 the quadratic potential 
stabilizes the split core boojum on large colloids, R/t, > 10, 
the single core boojum on small colloids, R/^ < 8, and the 
double core boojum at intermediate colloid sizes. By con- 
trast, in LCs with T] = —j, and the same quadratic potential, 
the split core is not observed and the single-double core transi- 
tion occurs at /?/^ ~ 45 (not shown). This behaviour is similar 
to the one described in previous sections, for quadratic poten- 
tials. When twist distortions are energetically unfavourable 
they are suppressed rendering the split core structure unstable. 
For a quartic potential with wi = W2 the split core is stable for 
both types of elastic anisotropy in the whole range of colloidal 
sizes. The size of the half-ring, given by the core splitting on 
the colloidal surface, increases with R and saturates at ~ 6£, 
for 77 = 2, and at < 10^ for rj ~ —j. Again, the quartic sur- 
face potential Eq. Q stabilizes split core boojums in LCs with 
negative elastic anisotropy. 

We end this section by noting that although, it may be dif- 
ficult to observe the reported boojum core structures using 
standard optical methods, their influence on the effective in- 
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FIG. 7: Defect core transitions driven by the colloidal size, for 
quadratic (circles) and quartic (diamonds) potentials with w ~ 37, 
at T ~ 0.16. For quadratic potentials, single core boojums are stable 
at small and spilt core boojums are stable at large colloidal sizes, R. 
The LC anisotropy T] = 2. For quartic potentials, the split core struc- 
ture is stable over the whole range of R. The solid lines are guides to 
the eye. 




FIG. 8: Schematic representation of two interacting boojum-coUoids. 
The far-field director is parallel to the z— axis. For large d pairs of 
boojums are aligned parallel to the far-field director. 



teraction between colloidal particles at short distances is sig- 
nificant, as will be discussed in subsection [231 The far-field 
nematic configuration, however, is not affected by the meso- 
scopic structure of the boojum cores. 



2.3. Two spherical particles 

In this section we shall calculate the two-body effective in- 
teraction potential f'^' {d, 9) ||45|] for colloids with a quadratic 
surface potential, Eq. (|6|l. For the definition of d and 0, 
see Fig. [8] f'^' is defined by decomposing the excess (over 
the free energy of the uniform nematic) free energy F2 = 
pLdG ^ fb (Qb)^ for the system of two particles as follows 

F2id,9)=2F'^^'^+F^-\d,e). (11) 
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FIG. 9: Effective interaction potential F^-\d) for several values of 
6. The surface potential is quadratic, Eq. (|6j, and the anchoring 
strength is w = vvj ~ 0.37. Elastic anisotropy rj =2, radius of the 
particles R ~ 6.6^, and reduced temperature T ~ 0.16. 

Here F^^^ is the excess free energy calculated independently 
for one isolated colloid. By definition F^^"^ tends to zero when 

In particular, we shall demonstrate how the non-linear ef- 
fects become dominant, leading to a qualitative change of the 
effective interaction, e.g., the potential F'^^ (d, 9=0) changes 
from repulsive at large d to strongly attractive at moderate to 
short distances. As we shall see this is driven by the complex 
behaviour of the interaction between the topological defects 
that results in a re-arrangement of their structure at short dis- 
tances. Analogous changes in the effective interaction have 
been re por ted experimentally (v^ and confirmed within LdG 
theory I120!] where it has been shown that the axial symmetry 
of the location of the boojum pairs is continuously broken, 
resulting in an effective colloidal attraction. 

Here we shall re-examine the effective colloidal interaction 
based on the results of the previous section where we have 
shown that the axial symmetry of the boojum cores may be 
broken at the level of a single colloidal particle. In what fol- 
lows we use the notation set in Fig.|8] The far- field director is 
aligned with the z-axis. Due to the spherical symmetry of the 
colloids, F*^' depends on the colloidal separation d and on the 
polar angle 9 only. In this section we consider colloidal par- 
ticles of equal size, R = 0.1 jim, the anisotropy of the elastic 
constants r] =2, the reduced temperature T ~ 0.16, and the 
other model parameters corresponding to the LdG parametri- 
sation of 5CB, described in subsection 12.11 For this set of 
parameters, the boojum cores exhibit a single core configura- 
tion over the whole range of anchoring strength w = w\. Ulti- 
mately, we aim at explaining the experimental results of Sma- 
lyukh et al. ITbIi . namely the force versus 9 curves, which dis- 
agrees strongly with the asymptotic quadrupole-quadrupole 
result, at moderate and short colloid-colloid distances. 

In Fig.|9]we plot F^^\d) in the weak anchoring regime. In 
this regime the interaction is repulsive when the inter-colloidal 
vector is either parallel or perpendicular to the far-field direc- 
tor, = or 7r/2 respectively, and is attractive at intermediate 
colloidal orientations, as expected for quadrupole-quadrupole 
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FIG. 10: Effective interaction potential F^^\d) for several values 
of 6. The anchoring strength is w ~ 37. Elastic anisotropy r\ = 2, 
radius of the particles R ~ 6.6t,, and reduced temperature T ~ 0.16. 
The insets depict the boojum configurations for several distances d 
and 6=0. In the attractive range, 2R <d < 3R, the inner boojums 
are ^ disclination rings represented by the iso-surface of the scalar 
order-parameter Q = 0.3Qi,. 



interactions. 

The character of the effective interaction changes signifi- 
cantly in the strong anchoring regime, as shown in Fig. [TO] 
For oblique, 9 = n/A-, and perpendicular, 9 = 7z/2, orienta- 
tions of the colloidal pair the interaction remains qualitatively 
the same: attractive when the colloids are at an oblique angle 
with the far-field director and purely repulsive when the col- 
loids are perpendicular to it. However, at = the effective 
interaction is no longer purely repulsive. At short distances, 
d < 3R, the particles start to attract each other This change 
results from a symmetry break of the boojum-pair configura- 
tion as shown in the insets of Fig. \TU\ At large distances the 
boojums are located at the poles of the particles, see inset (3) 
in Fig. [TOl and aligned parallel to the far-field director, as in 
the case of an isolated particle. As the separation between 
the colloids decreases, the repulsion of the inner defects in- 
creases. The change from repulsive to attractive colloidal in- 
teraction is driven by a change in the position of the inner 
defects, see inset (2) in Fig. [TOl This mechanism for attrac- 
tive colloidal interactions - the re-arrangement of the inner 
defects with quadrupolar symmetry - was reported almost 10 
years ago for 2D colloids with homeotropic anchoring ||8|. We 
note that although in the single particle configuration the boo- 
jums are in the single core state, the inner boojums undergo a 
transition to the split core state as the distance between the col- 
loidal particles decreases, see the inner boojums in Fig. [T0l 2). 
The symmetry break in the position of the inner boojums is 
accompanied by a symmetry break of the inner boojums core 
structure, which changes from the axially symmetric single 
core to the asymmetric split core configuration. The global 
minimum of F^^^ {d, 9) corresponds to close contact, d = 2R, 
at an oblique angle 9 close to 7r/4. This is in line with the ex- 
perimental results of iflslflTll where the boojum-colloids were 
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FIG. 1 1 : Effective interaction potential as a function of the distance 
d between two spherical colloidal particles, aligned with the far-field 
director, = 0. The surface anchoring potential is quadratic, Eq. 

Elastic anisotropy rj =2, radius of the particles R ~ 6.6^, and 
reduced temperature T ~ 0.16. The anchoring strength w varies in 
the interval [0.4,37]. The inset depicts the height of the free energy 
barrier as a function of w, exhibiting quadratic behaviour in the 
weak anchoring regime. 

reported to coalesce. In Sec. |3]we shall show that this is not 
the case for boojum-colloids in 2D systems. 

The effect of the anchoring strength w is trivial for colloidal 
orientations 6 = n/A^n/2, and thus we focus on = 0. In 
Fig.[TT]we plot F^"^^ {d, 0=0) for several values of anchoring 
strengths w. The dependence of the free energy barrier on the 
anchoring strength is illustrated in the inset. The barrier height 
is weakly discontinuous at w = w* ~ 0.7 delimiting two types 
of behaviour. For strong anchoring, w > w*, the free energy at 
distances close to the free energy maximum, ^f*, is a smooth 
function of d. This implies that as the colloidal separation de- 
creases the boojums rearrange smoothly, driving the change 
from repulsive to attractive interactions. For weak anchor- 
ing, w <w* , the nematic order at the poles of the particles 
is suppressed, but there are no defects. The system behaves 
almost linearly in this regime, i.e., Q 0= w and the resulting 
free-energy ^ w^. The effective interaction is repulsive at all 
d, as expected for two quadrupoles oriented at = 0. 

We conclude that the interaction between colloids aligned 
with the far-field director changes from repulsive to attractive 
at short-distances as the result of the re-arrangement of the 
positions and of the core configurations of the inner boojums. 

Finally, we shall give a detailed description of the angu- 
lar dependence of the effective force acting on the colloids. 
The experimental measurement of this force has been reported 
in Ref. ifisll . According to this reference brownian motion 
dominates the dynamics for d > 6R, but at shorter distances 
{d < 4-R) the dynamics is governed by the nematic mediated 
force. In this regime the force of the order of the elastic con- 
stant K has been found, ifisll . Since the measured elastic force 
was significant only for relatively small distances, strong devi- 
ations from the asymptotic quadrupolar interaction have been 
reported. For instance, the orientation corresponding to the 
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FIG. 12: Absolute value of the pairwise force as a function of 0, 
for several values of d. rj = 2, R ~ 6.6<^, T ~ 0.16 and anchoring 
strength w ~ 37. 



Strongest attraction has been found to depend on the distance 
between colloidal particles, while the interaction of quadrupo- 
lar type exhibits maximal attraction at « 49°, for any d. In 
addition, the attraction for 0" < < 70" and repulsion for 
75" < < 90" has been found, while quadrupoles repel at 
0=0. Finally, the measured radial dependence of the force 
(at fixed = 30") deviated systematically from the quadrupo- 
lar power law 0= d^^ for d < 3R. These findings suggest that 
non-linear effects dominate and that the superposition approx- 
imation does not hold in the range of distances that have been 
investigated. 

We have calculated the absolute value F of the elastic force, 
F = F^e, +Fgee, as well as its radial and polar components 
F, = —dFi^dcl dd and Fg = —{\/d)dF\^iQ/d6 by numerical 
differentiation. In Fig.[T2]the resulting magnitude of the force 
is plotted as a function of the relative orientation, 0, for sev- 
eral values of the particles separation d. We note that as the 
separation between the particles increases the magnitude of 
the force decreases and its structure as a function of de- 
creases. In Fig.[T2]the angular dependence of the force, at sep- 
aration d = 4R, resembles that of the quadrupolar force ll46ll 
with principal maxima at = and 6 = n. In addition, the 
force exhibits three secondary maxima at intermediate values 
of in line with the angular dependence of the quadrupole- 
quadrupole force ifisi l46ll . As the distance between the col- 
loids decreases the angular dependence of the force changes 
quite rapidly and quite drastically. Just below d ^3R the max- 
imum at = TT becomes a minimum (the force weakens) and 
the difference between the secondary maxima at intermediate 
values of becomes more pronounced. A similar angular de- 
pendence was observed in the force reported in Ref. lITsI at 
colloidal separations 3R ^ d ^ 4R. However, our theoretical 
results appear to be more sensitive to the inter-colloidal sep- 
aration and predict different angular behaviours at colloidal 
separations d = 3R and c/ = 4/? at = 0. At this orientation 
the experiments report a local minimum for colloidal separa- 
tions in the range 3R^d ^4Rf]J^. 

In an attempt to rationalise the discrepancies between the 
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FIG. 13: The radial component of the force between two colloidal 
particles as a function of 6, for several values of d. The values of the 
model parameters are the same as in Fig. 1121 
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FIG. 14: The polar component the force between two colloidal par- 
ticles as a function of the colloidal orientation, 6, at several values 
of the particles separation d. The values of the model parameters are 
the same as in Fig.ll2l 



theoretical and the experimental results we have plotted the 
radial and polar components of the force as a function of 9 in 
Figs. [T3]and[T4] The radial component of the force is found to 
be two orders of magnitude larger than the polar component, 
and thus it dominates the behaviour of the total force. This 
is in sharp contrast with the measured force, which is charac- 
terised by a radial component that differs significantly from 
the total force at short distances, around = (see Fig. 4 of 
reference ifisll . in particular the inset c) where the radial com- 
ponent of the force atd ~3R appears to vanish at = 0). The 
results for the polar component of the force plotted in Fig. [14] 
reveal that the angular dependence changes abruptly at = 0, 
by contrast to the behaviour of the force at = 7r/4 and 7z/2. 
This is related to the head-on interaction of the intermediate 
boojums and the associated discontinuity of Fg at = 0. This, 
probably, has had an impact on the experimental results as the 
colloidal configuration is unstable at = rendering an accu- 
rate force measurement, at this orientation, very difficult. Fi- 



nally, we note that the colloidal particles in the experiment are 
of the order of 2jj.m, while we considered colloids one order 
of magnitude smaller This may also account for some dis- 
crepancies between the theoretical and experimental results. 



3. TWO-DIMENSIONAL SYSTEMS 

Dispersed colloidal particles in smectic C (smC) films are 
realised either as inclusions of a lower order phase (isotropic, 
nematic, or cholesteric) or as smectic islands with a higher 
number of layers than the surroundings \2i 19} . Planar bound- 
ary conditions at the colloidal particle may lead to the nucle- 
ation of a pair of surface defects lfl9l l47ll . Planar anchoring 
may also lead to the nucleation of a satellite (bulk) defect ll44ll . 
which will not be considered here. At large distances the in- 
teraction between these particles is of quadrupolar type, and 
deviates from it at short-distances. By contrast to the obser- 
vations in 3D systems, where the spherical colloids with pla- 
nar anchoring come into contact and coalesce fll7il . inclusions 
in smC films maintain a well-defined separation deq — 2.1 R, 
where R is the inclusion radius. 

The distortions in smC films due to an isolated inclusion 
with planar anchoring, as well as the asymptotic effective 
pair interaction, have been described by invoking the elec- 
trostatic analogy and the superposition approximation. In 
the simplest approach lITsIl the distortions around an inclu- 
sion are modelled by the superposition of three solutions to 
the Laplace equation corresponding to one topological defect 
with charge -f 1 in the center of the inclusion and two — ^ sur- 
face defects. Then by using the superposition approximation 
the quadrupole-quadrupole effective interaction is obtained. 
However, this solution does not satisfy the boundary condi- 
tions at the inclusion boundary, and the superposition approxi- 
mation fails to describe the repulsion observed experimentally 
at small distances ifisll . A related approach ll2lll . based on the 
exact solution for the c— director (defined below) around an 
inclusion with arbitrary anchoring strength ll22l reveals that 
the repulsion appears simultaneously with the nucleation of 
the defects. Nevertheless, neither approach describes the non- 
linear effects responsible for the re-orientation of the defects 
at short distances. 

In what follows, we shall use the Landau description of 
smC films in order to calculate the interaction between two in- 
clusions with planar anchoring. We aim at a detailed descrip- 
tion of the defect structure, and we will show that the short 
distance colloidal interactions result from the re -arrangement 
of the defect positions, as in other 2D and 3D systems with 
quadrupolar symmetry iH . 



3.1. Free energy of smectic C films 

In smC phases the LC molecules are organized into 2D lay- 
ers with the local average molecular orientation n at a char- 
acteristic angle with respect to the layer normal v. When 
layer deformations are negligible there are no variations of 
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FIG. 15: Schematic representation of a smectic C. v is the normal to 
the layers, n is the average local molecular orientation, and c is the 
in-plane projection of n. 



potential, Eq. (|6]l, 

,,^,f«(Ll)\, ,16, 

JdQ. 2 V |Cz,| / 

where W is the anchoring strength, and V is the normal to the 
surface dQ. of the inclusion. We note that in 2D quartic terms 
do not change the surface potential qualitatively and will not 
be considered. Rescaling the variables we obtain 

f,.lf_«^(Ll)\l, (,7) 

J do. 2R V |c/,| / 

where we introduced the dimensionless anchoring strength 

w^WR/K. 



in the direction of v, V^^ = 0. In this case, the relevant order- 
parameter is the in-plane projection of the director n. This is 
an ordinary (variable length) 2D vector field c. 

The simplest mesoscopic model for smC films is obtained 
by expanding the free energy in terms of the invariants of the 
order-parameter c and its derivatives djCj ll48ll . The free en- 
ergy is then written as F = I jci{fb+ fe) where I is the 
thickness of the film, and the bulk and elastic free energy den- 
sities are, at lowest order. 



fb 



/, = |((V.c)2 + (Vxc)2 



(12) 
(13) 



where, for simplicity, we assume the one-elastic constant ap- 
proximation {Ki = ^3 = K). d{T) is assumed to be a lin- 
ear function of temperature T and is a constant. Dimen- 
sionless free energy densities (the free energy F has units of 
volume) are defined through the rescaling c — \Jb/ {2a)c and 
F = {b/a^)F, 



h = |C|2(|C|2-1) 



(14) 
(15) 



where f = \/bK/d{TY is the correlation length. Eq. (fl4l 
predicts non-zero values for the bulk scalar order-parameter 
\ct\ ~ l/\/2 at all temperatures. Topological defects corre- 
spond to regions of reduced orientational order |c| 0, and 
their topological charge is defined by the winding-number of 
the vector field c ll27ll . 

Due to the vector nature of the c director, it is possible to 
distinguish clockwise, anticlockwise or mixed planar anchor- 
ing on circular inclusions in smC films ||44[| . The orientation 
that is realised results from the interplay of the anchoring po- 
tential, the c director field in the vicinity of the inclusion, and 
the LC properties. We consider the case of mixed planar an- 
choring which favours the nucleation of boojum pairs 1I21I1 . 
the 2D counterpart of the 3D boojums investigated in subsec- 
tion 12.21 We use a surface potential of the Rapini-Papoular 
form [431 which is the 2D version of the quadratic surface 



3.2. Single circular particle 

In this section we consider a single circular inclusion in a 
smC film. The mixed planar anchoring is enforced by fix- 
ing the far-field c— director parallel to the y— axis. In this 
case a pair of boojums may nucleate at the particle surface. 
Within the FO elastic theory the c— director is a unit vector 
c = (cos<t>,sin<I>). The analytic solution for the orientational 
field, 4>, at arbitrary anchoring strength w, can be written as 



arctan 



[R/r) p{w) sin2^ 
1 - {R/r)~ p{w)cos,2(p 



(18) 



where r, 9 are polar coordinates, and p{w) = 
[2 / w) (^y/l+JwjTf- — \^ . The corresponding free en- 
ergy is Ffo = nlK{-\og{\-p^) + {w/2){\-p)). The 
solution, Eq. ( fTsT i, describes non-singular surface disclina- 

tions with a core size r^. = {2,R/w^)(^-s/\ /A - \ j iH, 

which for strong anchoring behaves as ~ AR/w. For 
any finite w the elastic free-energy is also regular, with the 
leading large w behaviour F{w ^ 1) ~ nKl\n{w). Singular 
topological defects appear only in the limit w — > 00. 

In the strong anchoring regime w — > 00 the exact solution, 
Eq. (fTsT l, may be approximated by the sum of a finite number 
of solutions <I>, to the Laplace equation, where <I>, represents 
a point singularity (or topological defect) of some winding 
number qi The number of these topological defects, their 
locations, and winding numbers, are chosen to comply with 
the boundary conditions. The behaviour of the c— director 
near the boojums suggests that these are — ^ topological de- 
fects. Since the colloidal particle with rigid planar anchoring 
accounts for a +1 virtual defect in the center, a naive guess 
would be to represent <J> as the sum of three point singulari- 
ties, which we write as (-5) + (+1) + (-5) 118]. Two (-5) 
singularities are located on opposite poles of the inclusion and 
are oriented along the far-field c— director However, a more 
careful analysis reveals that this ansatz violates the boundary 
condition at the surface of the inclusion, and additional virtual 
defects must be taken into account. The asymptotic w ^ °° 
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FIG. 16: Elastic free energy F''' — Fg, excess (over the free energy 
of the unifoitn system) free energy f''^, and surface free energy F, 
as a function of the planar anchoring strength w = WR/K. Points 
correspond to the numerical minimization of the Landau free -energy 
functional (Eqs. l ll4b . dlSb . QT}). Dashed lines are obtained using 
the exact solution, Eq. J18t . The radius of the particle R = 100<^. 



behaviour of the solution, Eq. (fTSl l. indicates that the correct 
ansatz corresponds to the triplet (— 1 ) + (+2) + (— 1 ), where a 
(+2) virtual defect is placed in the center of the inclusion and 
two (— 1 ) surface defects may be thought of as the sum of two 
(— j) defects, one of which is virtual and the other real llsill . 

In Fig. [T6]we compare the excess free energy F^^\ the sur- 
face free energy Fg, and the elastic free energy F*^'^ — Fs, ob- 
tained by numerical minimization of the Landau free energy 
functional (Eqs. (fl4l i. ( fTSl l, and iTT\ ) with the analytic results 
obtained by using the exact solution, Eq. (fTsT l, to the Frank- 
Oseen theory. Details of the numerical method may be found 
in the Appendix. In the weak anchoring regime the agreement 
is excellent. In the strong anchoring regime, however, the so- 
lutions start to deviate. Inspection of Fig. [T6]reveals that the 
discrepancies may be attributed to the differences in the sur- 
face free energies at large w. Indeed, the FO free energy is 
inaccurate in the presence of defects and ultimately diverges 
as ~ log(w). In the weak anchoring regime, w < 1, the excess 
free energy is dominated by the surface term which is linear 
in w, while the elastic free energy is quadratic, w^, as in 3D 
systems, see Fig.|4] 

Figure [17] illustrates the equilibrium c-director configu- 
rations for increasing values of the anchoring strength w. 
The grey-scale map corresponds to the optical transmit- 
tance, ///„ °^ sin^(2<t>), as seen in experiments under cross- 
polarisers. As the anchoring strength increases the inclusion- 
induced distortions, seen as bright fringes, grow. For w > 10 
a pair of topological defects nucleates at opposite poles of the 
particle (see FigfTTb) and d)). In order to quan tify the spatial 
extension of the boojums, we follow Ref. llSal and introduce 
an effective core radius as follows 
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FIG. 17: Equilibrium c— director configurations for anchoring 
strengths, a) w = 0.1, b) w = I, c) w = 10, and d)w = 100. The 
grey-scale map corresponds to the light intensity sin" 2<I> as seen 
in experiments under cross-polarisers. The radius of the particle is 
R 
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FIG. 18: Effective radius EQ- GU, of the 2D boojum core as a 
function of the anchoring strength w = WR/K. 



where the integral is over the region with one boojum. J as a 
function of the anchoring strength w is plotted in Fig. [18] We 
distinguish three regimes. As the anchoring increases the LC 
molecules align along the preferred surface orientation creat- 
ing small regions of low orientational order. At higher anchor- 
ing strength, I ^ w < 10, x °^ log(w), and the regions of low 
orientational order exhibit pre-nucleation of defects. Finally, 
in the nucleation regime, w > 10, the effective radius increases 
with w slowly, and eventually saturates at ~ 7/4^. 

By contrast to the 3D case, the 2D boojums are always 
point-like and never split. This is a direct consequence of 
the topological constraint which forbids half-integer defects in 
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smC films. In 3D nematics, however, half-integer ring discli- 
nations are not only topologically allowed, but are even ener- 
getically favourable ll34[ I35ll52l - l54ll . This is one of the rea- 
sons why single core boojums in 3D may split into half-ring 
disclinations. 



3.3. Two circular particles 

In this subsection we shall analyse the pairwise effective 
interaction between inclusions with planar anchoring in smC 
films. At d the inclusions "decorated" by boojums interact as 
2D quadrupoles H 

F,..MO)^'-^^P^. (20) 

This interaction is repulsive for — 7r/8 < Q — mK/2 < 7r/8 and 
attractive for 7r/8 < - mn/2 < 37r/8, with m = 0, 1, 2, .... 
At large distances the single-particle distribution of the c- 
director is only slightly perturbed by the presence of the other 
particle. As the separation decreases the perturbation in- 
creases and non-linear effects take over, changing the charac- 
ter of the interaction. The threshold distance d* at which the 
non-linear regime sets in depends strongly on the anchoring 
strength. 

We minimize the Landau free energy functional, Eqs. (fT4l i. 
(flSl l. ([TtI i numerically. Details may be found in the Appendix. 
In Figs.[T9]and|20]we plot the effective interaction free-energy 
F'^' = F + Fs — Fq — 2Fsi„^ig as a function of the particle sep- 
aration, d, for anchoring strengths w = 0.1, 1.0, 10, and 100. 
F +Fs — Fq contains the contributions to the free-energy from 
the distortions of the c— director, Fo = F{ci,) is the free en- 
ergy of the uniform bulk system, and F^ingu is the free energy 
of an isolated inclusion. At weak anchoring, w = 0.1, F'^' is 
always repulsive for perpendicular = 7r/2 or parallel 0=0 
orientations, and attractive for oblique orientations, = 7r/4. 

When the anchoring strength increases regions of reduced 
order develop around the poles of the particles leading to an 



FIG. 21: Equilibrium c— director configurations for interacting col- 
loidal particles with the anchoring strength w = 100 at relative orien- 
tation e = and separations a) d = 3.0R and b)d = 2.05R, R = 100|. 
c) close view of the configuration between the colloidal particles, for 
parameters as in b). The greyscale map coiTesponds to the light in- 
tensity oc sin^ 2$ as seen in experiments under cross-polarisers. 



increase of the elastic free energy, see Fig. [16] This also re- 
sults in an increase of the interaction strength at short dis- 
tances, as can be seen in Fig.[T9]for w = 1.0 and Fig. l20lfor 
w = 10. After the nucleation of topological defects the ef- 
fective interaction changes drastically as shown in Fig. |20] 
w = 100, at separations d < 4-R. In this regime the repul- 
sion becomes even stronger for orientations 0=0, n/2, but 
at = n/4 the free-energy develops a well defined local 
minimum at c/ ~ 2.4/?, which is followed by a repulsion for 
d < 2.4R. The latter prevents the inclusions from coalescing 
in agreement with the experimental observations lll9[|47l . We 
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recall that this is in sharp contrast with the behaviour of 3D 
boojum-colloids, where no such repulsion is observed. 

For strong anchoring and 0=0, F^'^\d) changes from re- 
pulsive to attractive, with a discontinuous slope, at ~ 2.45/?, 
seeFig.|20l w— 10. Similar behaviour is observed for w = 100 
at li ~ 2.10/?, and is not shown in Fig.[20]because the energy 
scale was chosen in order to emphasize the short range repul- 
sion and the local minimum for = k/A. Figure |2T| illus- 
trates the equilibrium c-director configurations for inclusions 
oriented along the far-filed c— director, = 0, and for the sep- 
arations a) d = 3.0/?, b) d = 2.05/?; the anchoring strength is 
w = 100. For d > 2.10/? {d > 2.45/? for w = 10) the inner 
boojums are aligned along the far-field c— director, giving rise 
to an overall repulsion between the inclusions. As the separa- 
tion decreases, the inner boojums slide away from the axis 
in opposite directions, as shown in Fig. I2T1:). This leads to 
a strong attraction between the particles. We emphasize that 
in 2D the re -arrangement of the inner boojums at = pro- 
ceeds discontinuously, while in 3D - continuously, see Fig.fTOl 
black circles. The contact value F^^^ {d = 2/?, 0=0) is, how- 
ever, larger than F^^' [d ~ 2.4/?, = 7r/4) which corresponds 
to the global minimum of F^^^ (c/, 0) for the case of strong an- 
choring. 



4. CONCLUSIONS 

We have carried out a detailed study of the effective col- 
loidal interactions for particles dispersed in a nematic host 
(3D) or in a smectic C film (2D). In both cases planar (de- 
generate in 3D) anchoring is imposed on the surfaces of the 
colloidal particles. Within the Frank-Oseen formalism in the 
strong anchoring regime, the boundary conditions are met by 
the creation of a pair of antipodal surface defects, boojums 
lfl6ll . For an isolated particle the vector connecting two boo- 
jums aligns itself with the far-field director, as required by 
the global boundary conditions UtIi . The resulting director 
field has quadrupolar symmetry, and the ensuing effective col- 
loidal interactions exhibit quadrupole-quadrupole asymptotic 
behaviour (dHIiHIl, both in 2D and in 3D. However, at 
short distances, where the superposition approximation fails 
and the non-linear effects dominate, the experimental results 
report significant deviations from this quadrupolar behaviour 
The most important of these is the crossover from a repulsive 
to an attractive interaction at some threshold distance for parti- 
cles aligned along the far-field director, = 0. The threshold 
distance is ~ 3/? in 3D, see Fig. [TOl and ~ 2.5/? in 2D, see 
Fig. |20] This crossover is driven by the configurational re- 
orientation of the inner boojums as the distance between the 
particles decreases. 

We have shown that the short-range colloidal interactions 
in 3D result not only from the re-arrangement of the defect 
positions, as predicted in 2D a decade ago fl] and described 
recently in 3D ll20ll . but also from changes in the structure of 
the boojum cores. The description of the structural transitions 
between boojum cores is a challenging theoretical problem as 
distinct structures, characterized by tensor order-parameters 
that vary in regions of the order of the bulk correlation length. 



have similar free energies. Their stability results from a deli- 
cate balance of various contributions and structural transitions 
may be driven by the anchoring strength, the temperature or 
the colloid radius. We have concentrated on the mesoscale and 
addressed the structure and dynamics of the boojums and the 
resulting colloidal interactions based on the LdG free energy 
in 3D and 2D, for spherically symmetric colloids. We have 
used finite elements methods with adaptive meshes in order to 
resolve the structure of the defect cores and establish the na- 
ture of the short-distance effective colloidal interactions. We 
have found that the defects are point like in 2D, but acquire a 
rather complex structure in 3D, which depends on the combi- 
nation of the anchoring potential, the colloid radius, the tem- 
perature and the LC elastic anisotropy. In addition, we have 
found defect-core transitions between ;) single core, //) double 
core, and ;//) split core structures of the boojum pairs. 

Single core boojums are uniaxial along the z— axis (see 
Figs.[T]5!,|2]5!). The scalar order-parameter vanishes at one point 
above the colloidal surface, where the system is isotropic. Be- 
low this point the scalar order-parameter is negative and the 
nematic order is oblate. Above that point the scalar order- 
parameter is positive indicating that the nematic order is pro- 
late. By contrast, double core boojums are biaxial along the 
z— axis (see Figs. [T}7, [JJ?). The surface of maximal biaxiality 
intersects the z— axis at three distinct points separating regions 
of low nematic oblate and prolate order Split core boojums 
are also biaxial along the z— axis. The surface of maximal bi- 
axiality intersects the z— axis at two distinct points that delimit 
the region of oblate nematic order. Figs. [1];, [2j;. As a gen- 
eral conclusion we find that on small colloids boojum cores 
are axially symmetric (see Fig. |7), point-like with index +1 
- the index referring to the charge of the projected director 
field on the surface -, which are stable at high temperatures 
(see Fig. |6]l and relatively weak anchoring (see Fig. |5). On 
large colloids at low temperatures and strong anchoring, the 
axial symmetry is broken and the boojum +1 point-like cores 
split into a pair of +j point-like surface defects, connected 
by a disinclination line (see Fig. [T]:). A structure without a 
fully developed disinclination line, the double core boojum, 
was also found (see Figs.[Tl7,|3]l. The detailed structure of the 
cores as well as the transitions between the different configu- 
rations depend in detail on the colloid and LC parameters. We 
stress that the details of the core structure do not affect the far 
field configuration, i.e., on a spherical surface surrounding the 
boojum cores the order-parameter distribution resembles that 
of a point boojum. Their effect, however, on the short-distance 
interaction may be significant. Indeed, we have shown that the 
short-distance attraction at = 0, which results from the re- 
orientation of the inner boojums, as the distance between the 
colloids decreases, corresponds to a structural change from 
single core to split core for the inner boojums pair. 

We have shown how the non-linear effects become domi- 
nant, both in 2D and 3D, leading to a qualitative change of the 
effective colloidal interaction, namely, the crossover from the 
large-distance quadrupole-quadrupole repulsion for a range 
of orientations around = 0, to an attraction at moderate to 
small distances, see Figs.[T0l|20] This results from the com- 
plex dynamics of the topological defects that rearrange their 
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position at short distances both in 2D and 3D. This attraction, 
along the far-field direction, agrees with experimental obser- 
vations H [H [H, m, and with earUer calculations S 
where it is shown that the axial symmetry of the boojum-pairs 
is broken when the repulsion changes into the attraction. The 
mechanism for attractive colloidal interactions - driven by the 
re-arrangement of defects with quadrupolar symmetry - was 
re|)orted in 2D for homeotropic anchoring nearly a decade ago 

We emphasize an important difference between the 2D and 
3D effective pair interaction potentials. In 3D the global min- 
imum occurs at contact d ^2R for orientations at oblique an- 
gle with the far-field director In 2D the optimal orientation of 
two inclusions is also oblique, but there is a free energy barrier 
~ 0.51K keeping the particles apart. This distinction is in line 
with experimental observations in 2D and 3D llT3[fl8[fl9[l23ll . 

Finally, we have calculated the force between boojum col- 
loids in 3D, see Fig. [12] We have found, that as the dis- 
tance decreases the angular dependence of the force changes 
rapidly and quite drastically. Just below d = 3R the maximum 
at = TT becomes a minimum and the difference between 
the secondary maxima at intermediate values of 6 becomes 
more pronounced. This angular dependence was actually ob- 
served by ioil at larger colloidal separations, 3R < d < 4-R. 
The theoretical results appear to be too sensitive to the inter- 
colloidal distance, when compared to the experimentally mea- 
sured forces. 
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Appendix A: Numerics 

In this section we discuss briefly how the scheme for nu- 
merical minimization of the free energy functionals is imple- 
mented. In order to discretize the continuum models we resort 
to a finite element method llssll with adaptive meshes. For 3D 
systems the nematic and the colloidal particles of radius R are 
confined in a cubic box of linear size 30 x R. For smectic C 
films the system is confined in a 10/? x 10^ square box. 

1. Meshing 

The surfaces of the spherical colloidal particles are dis- 
cretized using the open source GNU Triangulated Surface Li- 
brary (5^. The library uses a recursive subdivision algorithm 
in order to triangulate the surface of a unit sphere. Starting 
from the icosahedron as the first approximation the next lev- 
els of refinement are constructed by subdividing each triangle 



into four smaller triangles. We have used the sixth refinement 
level corresponding to a surface mesh with 20480 triangles. 

In the next step the domain Q. accessible to the nematic is 
triangulated. For 3D systems the triangulation is carried out 
using the Quality Tetrahedral Mesh Generator, TetGen, iIstIi 
library. TetGen generates boundary conforming Delaunay 
llssll meshes of 3D domains with a piecewise linear bound- 
ary. The boundary can also be specified as a surface mesh, 
which is the case for the triangulated surfaces of the spherical 
particles. TetGen supports an isotropic size conforming trian- 
gulation, i.e., the resulting elements conform to a given local 
size map. For 2D systems the corresponding nematic domains 
are triangulated using the INRIA's BL2D software |59fl. The 
BL2D package supports the creation of anisotropic meshes, 
where the triangulation process is governed by specifying a 
metric map. 



2. Minimization 

We use linear elements both in 2D and 3D, i.e. the values 
of the nematic field Q, are specified only at the vertices x, of 
the mesh and a linear interpolation is used in order to deter- 
mine the values of Q(x) at some other point x of Q.. Due to 
the discretization the volume and surface integrals in 3D are 
replaced by sums of integrals over tetrahedral and triangular 
elements, respectively, and in 2D by sums of integrals over 
triangular and line elements (the discretized perimeters of the 
colloids). Each integral is evaluated numerically using gen- 
eralized Gaussian quadrature rules for multiple integrals 1I6OI1 . 
For integration over tetrahedra a fully symmetric cubature rule 
with 1 1 points, which is exact for polynomials of degree 4, is 
used [60^.61,1 . For integration over triangles a fully symmetric 
quadrature rule with 7 points, which is exact for polynomials 
of degree 5 is used. Finally, one-dimensional integrals 
are evaluated using the 5 points Gauss-Legendre quadrature, 
which is exact for polynomials of degree 9. 

During the minimization the values of the tensor order- 
parameter Q at the vertices of the bounding box are kept fixed, 
equal to the values of the uniaxial bulk nematic. The values of 
Q; at all other vertices are obtained by numerical minimiza- 
tion. We used the INRIA's M1QN3 |63] optimization routine. 
The routine implements a limited memory quasi-Newton tech- 
nique (limited memory Broyden-Fletcher-Goldfarb-Shannon 
method) of Nocedal Q). 



3. Adaptive mesh refinement 

The systems are characterised by two widely different 
length scales. One is given by the nematic correlation length 
and the other by the radius of colloidal particles. Therefore, 
we resort to adaptive re-meshing technique in order to obtain 
a sufficiently good approximation for the "exact" solution in 
a reasonable amount of time and with the finite computing re- 
sources at hand. The main objective of the adaptive refinement 
is to construct a mesh which meets the given interpolation er- 
ror tolerance for a minimal number of elements. Such meshes 
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are called optimally-efficient JHH] . Let / be one of the com- 
ponents of the exact solution Qij to our problem, and /' its 
linear interpolation on some mesh. It can be shown that the 
interpolation error associated with some element E/t (either 
a tetrahedron or a triangle) satisfies the following inequality 
M 

max|/(x)-/'(x)| ^ ma\\xaJ^aRi'>^o)xR\, (Al) 

where xq is some vertex of the element Ei^, x is the position 
within element Ei^ measured relative to xq, and J^ap (^o) = 

gj^^g^^ |xo is the Hessian of the exact solution / evaluated at 
the vertex xq. If we define the absolute value of the 2D Hes- 
sian as 

W-'^'Co'' |°,|)^. 

where is an orthogonal matrix that diagonalizes J^, then 
the estimate for the interpolation error on the element £4 can 
be written in a simpler form 

max|/(x)-/'(x)| </2fmax(Ai,A2)- 

hii is the diameter (the length of the longest edge) of the ele- 
ment Ell. It is straightforward to write down similar expres- 
sions for the 3D case. Therefore, if the sizes of the elements 
are chosen such that /!|max(Ai , A2) — const for all the ele- 
ments, then the resulting interpolation error will be approxi- 
mately equally distributed among all the elements. 



As it was mentioned above the BL2D package supports 
anisotropic meshes. The triangulation is governed by specify- 
ing a metric map and the resulting mesh is characterised 
by edges of unit length according to the metric ^ ll66ll . We 
assume that the required metric is proportional to the Hessian, 
^ 0= 1^1- In order to estimate at some vertex xq) we use 
the weak definition of the Hessian ll66ll 

where the integral is over the elements which share the vertex 
xo, and 0" is the piecewise linear hat-function associated with 
vertex 0: 0" = 1 at xq and 0" = at any other vertex. We use 
linear interpolation in order to obtaine at some internal 
point X of an element. 

TetGen implements an isotropic mesh refinement strategy 
which is based on a maximum local element volume con- 
straint. The constraints on the volume of elements are ob- 
tained by applying the equidistributing principle, where the 
nodes of the mesh are chosen such that for each element Ej^ 
the following condition holds: 

/ ^/\deU^\d''r = const. (A2) 

In Ref. ll65ll . it is shown that an optimally-efficient 2D triangu- 
lation which minimizes the interpolation error fulfills asymp- 
totically the equidistributing principle (IA2b . We assume that 
the same is valid in the case of 3D triangulation. 
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